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Abstract 

We present models of dormancy in a planktonic culture and in biofilm, 
and examine the relative advantage of short dormancy versus long dor- 
mancy times in each case. Simulations and analyses indicate that in plank- 
tonic batch cultures and in chemostats, live biomass is maximized by the 
fastest possible exit from dormancy. The lower limit of time to reawaken- 
ing is thus perhaps governed by physiological, biochemical or other con- 
straints within the cells. In biofilm we see that the slower waker has a 
defensive advantage over the fast waker due to a larger amount of dormant 
biomass, without an appreciable difference in total live biomass. Thus it 
would seem that typical laboratory culture conditions can be unrepre- 
sentative of the natural state. We discuss the computational methods 
developed for this work. 

1 Introduction 

Microbial populations, particularly those in biofilms (sessile, matrix encased 
communities, see [IT] for an overview), can contain cells in varying phenotypic 
states. An important difference between planktonic (free-swimming) and biofilm 
environments is that the former is generally well-mixed whereas the latter is un- 
mixed and spatially heterogeneous as a result. As a result of their self-generated 
spatially variable environmental, biofilms demonstrate spatially diverse ecology 
P01 [TT] . Such diversification may be advantageous for defense against an un- 
certain and temporally varying environment. For example, though cell states 
that are more tolerant to antimicrobial challenge may be less competitive in the 
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absence of that challenge, their presence can improve community survivability 
against attack. 

Here we consider one such defense mechanism, dormancy (possibly related to 
the phenomenon of persister cells [SJ[T2]) where, in response to an environmental 
stress, cells differentiate into a protected, slow- or non-growing condition [71 [T3]. 
Bacteria in planktonic states have been found to revive faster from dormancy 
than those in a biofilm state [7J- Thus it would seem that dormancy- regulating 
parameters are subject to influence of environmental variation, at least of the 
sort found in biofilms. Here we wish to use modeling tools in order to gain insight 
into role and regulation of dormancy in spatially mixed systems (batch and 
chemostat microbial communities) and unmixed systems (biofilm communities). 
Our attention is directed to the relative advantage of short dormancy versus long 
dormancy times in the cases of batch, chemostat, and biofilm states. 

In the process, we also present computational tools designed to study dor- 
mancy within batch, chemostat, and biofilm population dynamics, in particular 
with respect to competitiveness. These tools are an extension of those discussed 
in [5] for persistence and senescence, primarily in the numerical methods used 
to solve the more general physiological structure used in this paper, and we ex- 
pect they will have wider applicability to descriptions of physiological states in 
both mixed and unmixed microbial communities . The physiological structure 
is represented by a continuous variable. Compartmentalized dormancy models 
have been considered elsewhere [131 [Tl] . 

This paper is organized as follows. We present models of chemostat and 
batch cultures, and asymptotic analyses of their long-time behavior. We then 
derive the biofilm model. We compute numerical solutions of the model equa- 
tions for the batch, chemostat and biofilm cultures, and discuss the numerical 
methods developed for these computations. We conclude with the implications 
of our results. 

2 Models of Dormancy in Chemostat and Batch 
Cultures 

We introduce s € [sq, s*] to index the dormancy state of individual cells, with 
So the value at which cells enter dormancy and s* the value at which cells leave 
dormancy and become active. Cells progress through dormancy states with 
"speed" g(s,c), where c is concentration of relevant chemicals (e.g. substrates 
or antimicrobials); for example, large concentrations of substrates and/or small 
concentrations of antimicrobials imply larger value of g. While dormant, cells 
do not grow and divide; on the other hand, dormant cells are presumed to be 
hardier in response to environmental stress. 

Let t > represent time. Let u(t) represent the density of active cells, v(s, t) 
represent the density of dormant cells, and c(t) be a vector of substrate chemical 
species concentrations. Let the operator d y denote partial differentiation in the 
subscript variable y. The active cell population is modeled by an ordinary 
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differential equation, t > 0, 



^juit) = b(c)u(t) - fj, u (c)u(t) 

cell division death 

+ g(s*,c)v(s*,t) - h(c)u(t) -d u(t), (la) 

exit from dormancy entrance to dormancy washout 

where do is the chemostat dilution rate (roughly, inverse of the time scale for 
the chemostat contents to be flushed). The functions b and (i u account for cell 
division and cell "death" . We use "death" as shorthand for all forms of inertness 
not tied to strategic dormancy on the part of the bacteria. The function h > 
is the dormancy rate of active cells. Let g > (as above) and /i v > denote the 
reactivation and death rates, resp., of the dormant cells. We use a physiologically 
structured equation for the dormant cell population, 

d t v{s, t) + d s (g(s, c)v(s, f )) = -fi v (s,c)v(s, f) - d v(s, t) , (lb) 

reactivation kinetics death washout 

g(s ,c)v(s ,t) = h(c)u(t) , (lc) 

creation of newly dormant cells 

for so < s < s*. For the substrate chemicals, we have 

^c(t) = - f(c(t),u(t), V(; t)) + d (C (t) - c(t)), (Id) 
at ^ v ' v v ' 

chemical reactions chemostat turnover 

where / is the vector of reactions and Co(t) is the input concentration vector 
from the chemostat tank. Initial conditions are it(0) = Uq, v(0,s) = fo(s), an d 
c(0) = c . The model for a batch culture is obtained from (JXJ) by setting do = 0. 



3 Long-Time Behavior 

In this section we examine the long-time behavior of chemostat models for steady 
and periodic cases. 

3.1 Steady Chemostat 

We consider the long-time behavior of the steady (Co(t) = Co) chemostat system 
by studying the time-independent solution of equations (JTJ) . In the steady state 
c(t) = c, so that we can define a new dormancy coordinate a G [0, o*] by 

f S ds' 

with a* = a(s*). Setting the time derivative to zero, (jlb[) together with (flcl) 
can be solved to obtain 

„( ) = J^}]L e -doa e -f ^ v (a',c)da' (3) 

v ' 9(a*,c) v ' 
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where u is the steady state value of the active cell density and s has been 
replaced by a. The first exponential factor accounts for loss due to washout 
and the second for loss due to death. Plugging into (fTa)) and again setting the 
time-derivative to zero, we obtain 



If u ^ 0, then the second factor of (QJ provides an equation for c. Assuming 
that fi u , fi v , and h are all decreasing functions of c, we can write that second 
factor in the form 6(c) — a — a(c), where a is a constant and a(c) is a decreasing 
function of c such that a(c) — > as c — !• oo. Thus, assuming that 6(c) is a 
monotone increasing function of c and that 6(c) > a for c sufficiently large, then 
it follows that there is a unique value of c that solves (Q| with u^O. However 
if that value is larger than Co or if 6(c) < a for all c, then the only admissible 
solution of (@]) is u = (washout). 

Finally, given the solution for c, equation (|ld[) can be solved to obtain the 
long-time behavior of u by setting the time derivative to zercP|. In the case that 
/ is monotone increasing in its arguments, we obtain a unique solution for it. 

For two (or more) species competing in the same chemostat, the one that has 
the steady-state solution with smallest value of substrate c is the only long-time 
survivor |16j ; it excludes the other species by continually reducing substrate until 
substrate level is below the others' steady-state requirements. Thus, from (j4j it 
is apparent that the smaller the size of h (i.e., the lesser the likelihood of entry 
to dormancy) and the smaller the value of a* (i.e., the shorter the dormancy 
period), the more competitive the species. Note that a species which does not 
go dormant at all will outcompete an otherwise similar species which does. 

3.2 Periodic Chemostat 

We consider next the long-time behavior of a periodic chemostat with input 
substrate concentration Co(ujt) where Co(a + l) = Co (a), and offer asymptotics 
for two special cases. 

3.2.1 Fast Oscillations 

We suppose that the chemostat oscillation period is short compared to all other 
time scales of interest. In this limit, the chemostat will oscillate many times 
before the microbial inhabitants can react. This intuition suggests a multiple 
time scale expansion with a slow time t\ and a fast time ti defined by 



1 I{ u = 0, then JTdJl requires that do(Co — c) = /(c, 0, 0); generally /(c, 0, 0) = in which 
case c = Cq . 




k = t, 



t 2 =e-H, 



(5) 
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where e = ui 1 <C r for any inherent time scale r in the system. Note Co = 
Co(i2)- We expand 

u = u (h,t 2 ) + «ii(ii,t 2 ) + . . . , (6a) 

v = v (s,t 1 ,t 2 ) + evi{s,ti,t 2 ) + ... , (6b) 

c = c (t 1 ,t 2 ) + eci(t 1 ,t 2 ) + . . . , (6c) 

We suppose that the solution (u,d,c) approaches periodicity with period 
oj^ 1 for long times, so we look for a solution independent of slow time t±, i.e., 

u = u (t 2 ) + eui{t 2 ) + . . . , (7a) 
v = v (s,t 2 ) + evi(s,t 2 ) + ■ ■ • , (7b) 
c = c Q {t 2 ) + ecx{t 2 ) + ... , (7c) 

In this case, d/dt = e~ 1 d/dt 2 . Then to lowest order (= 0(e -1 )), system |T]) 
becomes 

duo _ dv^ _ 3cq _ 
dt 2 dt 2 dt 2 
so that uq, cq are constants and vq = vq(s), i.e., 



(8) 



u = uq + eui(t 2 ) + . . . , (9a) 
v = vo(s) + evi(s,t 2 ) + ■ ■ ■ , (9b) 
c = c + eci(t 2 ) + . . . , (9c) 

At the next order (= O(e )), system (JT|) becomes 

-^—ui{t 2 ) = b(c Q )u - fi u (c )u + g(s* ,c )v (s*) - h(c )u - d u , (10a) 
dt 2 

d 

-^-vi(s,t 2 ) = -d s (g(s,c )vo(sj) - fi v (s, cq)vq(s) -d v (s), (10b) 
ot 2 

d 

-^—cx{ti) = -f(co,u ,v (-)) + do(Co(t 2 ) - c ), (10c) 
ot 2 

with g(so, co)w(so) = h(co)uo- Averaging over a period u , we obtain 

= b(c a )u - fx u (c )u + g(s*,c )v (s*) - h(c )u Q - d u , (11) 
= -d s (g(s,co)v (s)) - fi v (s,c )v (s) - d v Q (s), (12) 
= -f(co,ua,v (-)) + d (C - c ), (13) 

with g(so, cq)v(sq) = h{co)uo, where Cq is the average of Co(t 2 ) over one chemo- 
stat oscillation period. This system, the same as was solved previously in the 
steady chemostat case except with Cq replacing Co, has essentially the same so- 
lutions for Mo, fo(s), and Co; the next order terms ui(t 2 ), Vi(s,t 2 ), c\{t 2 ) add a 
correction of O(e). Note thus that the same conclusion holds: a species without 
dormancy will outcompete an otherwise similar species which can go dormant 
(because the fast oscillating chemostat acts like a steady chemostat with input 
substrate Co). 
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3.2.2 Slow Oscillations 



We suppose now that the chemostat oscillation period is long compared to all 
other time scales of interest, i.e., that the chemostat can nearly reach equilibrium 
before input Co(wt) changes noticeably. Intuition again suggests a multiple time 
scale expansion with a slow time t\ and a fast time t 2 , in this case defined by 

h = et, t 2 = t, (14) 

where e = ui <C t _1 for any inherent time scale r in the system. Note that 
Co = Cb(ti). We expand 

u = u (ti,t 2 ) + eui(ti,t 2 ) + . . . , (15a) 
v = v (s,ti,t 2 ) + evi{s,ti,t 2 ) + ... , (15b) 
c = c {ti,t 2 ) + eci(t 1 ,t 2 ) + ... , (15c) 

We suppose quasi-equilibrium in the sense that Uq, vq, and cq are indepen- 
dent of fast time t 2 . Noting that d/dt = ed/dti + d/dt 2 , then at its slowest, e° 
order, system ([T]) becomes 

= b(c )u a - li u (c )uq + g(s*, cq)v (s*) - h(co)u - d u 0l (16a) 
= -d s (g(s,co)v ) - fj, v (s,c )v - d v , (16b) 
= -f{co,u ,v (-)) + d (C {ti) - c ), (16c) 

with g(so, co)v(sq) = h(co)uo. Note that t\ is essentially a parameter, appearing 
explicitly only in the input substrate concentration Co(ti). Thus to zeroth order, 
the quantities u, v, and c obey system (|16[) . which is the same as the steady 
chemostat except with parameterized input substrate. Hence we again conclude 
that a species that does not go dormant will outcompete an otherwise similar 
one that does. Note one caveat though: if G{t\) drops below the minimum 
required to sustain a particular population (see Section I3.1[) at any point in its 
cycle, then extinction may occur. 

4 A Model of Dormancy in a Biofilm 

For the biofilm model we remove the chemostat-specific terms and extend the 
system ([1]) to include a spatial domain SI consisting of stratified subdomains 
B t for biomass and tt\B t for the bulk fluid. There are two moving interfaces 
in SI: r t separating B t from the rest of fl, and a bulk-substrate interface 
that is a fixed height Ht above T t . The biofilm rests on a surface, denoted by a 
lower boundary, Tb- The spatial domains are illustrated in Figure [TJ The active 
and dormant cell populations, and the chemical concentrations, now depend on 
x e SI. Conservation of biomass yields, for t > and so < s < s* , 

d t u(x, t) + V • J M = 

6(c)m(x, t) — jjL u (c)u(yL, t) + g(s*, c)w(x, s* , i) — /i(c)u(x, i), (17a) 
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Figure 1: Spatial domains for the biofilm model. 

d t v(x,s,t) + d s (g(s,c)v(x,s,t)) + V • J v = -fi v (s,c)v(x,s,t), (17b) 

g(s , c)v(x, s , t) = h(c)u(x, t), (17c) 

d t w(x,t) + V • J w = /x u (c)u(x,i) + / fjb v (s,c)v(x,s,i) ds, (17d) 

J SO 

9 t c(x, t) + V • J c = /(c(i), u(x, i), v(x, •, i)), (17e) 

(17f) 

where V- denotes divergence in space, 3 y denotes the flux of subscript variable 
y, and where we assume appropriate initial conditions and boundary condition 
on the spatial domains. 

Assuming Fick's Law gives J c = —DVc for diffusion constant D. The sub- 
strate masses are also subject to advection, but the velocity is sufficiently slow 
that we can neglect the advective contribution to the flux. Likewise, substrate 
material diffusion time scales are at least several orders of magnitude larger 
than the those at which bacteria grow or advect, allowing us to make a quasi- 
steady-state assumption so that 

-DV 2 c = f. (18) 

Let i?(x, s, t) and p(x, s,t) denote the volume fraction per dormancy state 
and density per dormancy state relative to volume fraction, resp., of dormant 
cells. We assume incomprcssibility of biomass with p(x,s,t) = p* for positive 
constant p*. We also assume, based on the fact that the main constituent of 
all cells is water, that active and inert cells have the same incomprcssibility 
properties, and the same densities relative to volume fractions, p* , as dormant 
cells. We let z/(x, t) and r)(x,t) denote the volume fractions of active and dead 
cells, resp., which are related to the density of active and dead cells by u — p*v 
and w = p*r/. 

Assume the biofilm polymer matrix exists in proportion to cell density. We 
require the biomass volume fractions to total to one so that 

z/(x,t) +r?(x,f) + [ i?(x,s,t) ds = 1. (19) 

J So 
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Assuming that transport of biomass, including dormant cells, is governed 
by an advective process, with a volumetric flow u(x, t) for all classes and ages, 
gives the fluxes J u = p*uu, J v = p*vu and J w = p*wu. As in [5], we follow 
[TJ [TU] and assume that the volumetric flow is stress driven according to 

u = -AVp, (20) 

where p(t, x) is the pressure, A > the Darcy constant, and p — in il\B t . Pres- 
sure is determined in order to enforce incompressibility in response to growth 
and hence (|2"0|) can be viewed as a balance of growth-induced stress against 
friction. Other choices of force balance are possible. 

Substituting u — p*v, J u = p*uu, v — p*d, J v = p*vu, w = p*rj, and 
Jw = p*wu into equations (fTTf gives 

d t v{x,t) +V - (up) = 

b(c)v{x,t) - fi v (c)u(x,t) + g(s*,c)&(x, s*,t) - h(c)v(x,t), (21a) 

9 t i?(x, s, t) + d s (g(s, c)0(x, s, t)) + V • (u0) = -/x„(s, c)tf(x, s, t), (21b) 

g(s , c)i?(x, s , i) = h(c)v(x, t), (21c) 

d t r](x,t) + V • (u?y) = p u (c)v(x,t) + / p v (s,c)tf(x,s,t) ds, (2 Id) 



so 



-DV 2 c = /, (21e) 



with appropriate initial conditions and boundary condition on the spatial do- 
mains. 

Integrating (|21bj) over s gives 



#(x,s,i) ds +g(s*,c)i?(x, s*,t) -p(s ,c)i?(x,s ,t) 

so / 



— hu 



-dtv—dtV 



+ V (u/ ti(x,s,t)d S )=- [ p v (s,c)'d(x,s,t) ds. (22) 

\ •'So / •'so 

=V-u(1-j/-J7) 

Substituting for — dtv — dtf] yields 

V ■ u = b{c)v. (23) 
Substituting u = — AVp gives an equation for the pressure in Bt, 

- AV 2 p = b(c)v. (24) 
Distributing the divergence operator gives 

V • (uv) = —AVp • Vv + b(c)v 2 , (25a) 

V • {ad) = -AVp • V?9 + b(c)v8, (25b) 

V • (urj) = -AVp • V-q + b(c)vq. (25c) 
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We see from ([Ml) that p is proportional to A -1 , so that XV p is independent 
of A. Consequently v, fl, and rj are independent of A, allowing us to set A = 1. 

We impose periodic and other boundary conditions, similar to what was done 
in pQ , to obtain the complete model. The active cell volume fractions satisfy 

<V(x,t) ~Vp-Vv = -fi v (c)u(x,t) + #(s*,e)tf(x, s*,t) 

- h(c)u(x, t) + b(c)i/(x, t)(l- 2/(x, *)) , (26a) 

for x G B tl t > with conditions 

-^ = o, xer B ,i>o, (26b) 

oz 

Kx,0) = i/ Q (x), xeB ( , (26c) 

where z denotes the spatial variable orthogonal to the surface Tb , and vq is the 
initial active cell population. The dormant cell volume fractions satisfy 

S t tf (x, a, t) + d s (g(s, c)0(x, s, *)) - Vp • W = 

- n v (s,c)i}(-x.,s,t) - 6(c)i/(x,i)i?(x,s,t), (26d) 

for x £ B t , s > So,t > 0, with conditions 

c?(s , c)i?(x, s , *) = fe(c)i/(x, t), xeB(,i>0, (26e) 

— = 0, x G T B , t > 0, s > s , (26f) 
oz 

tf(x,s,0) = 0, xeSt,s>s . (26g) 
The fully inert cell volume fractions, including necrotic cells, satisfy 

d t ri(x,t) - Vp • Vrj = 

/tiu(c)i/(x,t) + / p v (s,c)$(x,s,t) ds - 6(c)i/(x,t)i7(x,t), (26h) 

for i £ Bj, t > 0, with conditions 

^ = o, xer B ,f>o, (26i) 

? 7 (x ) 0)=77o(x) ) x£B () (26j) 

where 770 is the initial inert cell population. Pressure satisfies 

-V 2 p = b(c)v, xeB t ,t>0, (26k) 
p = Q, xer t ,i>0, (261) 

^=0, xer B ,(>0. (26m) 
Oz 
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Let / = [fx , . . . , f m ] and c = [c± , . . . , c m ] . The chemical species satisfy, for 
j = l,...,m, 

-DjV 2 ^ =fj, x G 0,* > 0, (26n) 
£ - 0, x G fi\B t , (26o) 

cj = c*, xer fltI f>o, (26 P ) 

-k L =Q, x£r B ,i>0, (26q) 

where the Z)j are chemical diffusion coefficients and the c* are the chemical 
concentrations in the bulk fluid. The normal velocity of the interface is 
given by 

-y p -n=-^ > (26r) 
where n is the unit outward normal of Tg . 



5 Computations 

In this section we present computational results for models of batch cultures, 
chemostat cultures, and biofilms. 



5.1 Batch Culture Dormancy As a Response to Nutrient 
Depravation 

We let c(t) be a scalar value representing nutrient. We choose the functional 
forms 

6(c) = (27a) 
h(c) = ^-+e h , (27b) 

5(S)C)= V (27c) 
7 + c 

M«(c) = M«, (27d) 

Hv(s,c)=0, (27e) 



f(c, u, v) = -y^—j (u + f (e- s ^ + e-C 1 -)/*.)* ds\ 



(27f) 



«o(«) = 0, (27g) 

with constant fi u , rate constants k, kh, k g , saturation constant 7, Monod con- 
stant (, and yield constant Y. We take the baseline parameters k = l/4hr, 
7 = 4g CO DB/ m3 > eft = 0-05, M« = 0.005/hr, and Y = 0.63g CODB /g CODS [19]. 
The units gcoDB an d Scods are the chemical oxygen demand of biomass and 
substrate mass, resp. We assume a small £, say £ = 7/2O = g CO DB/5 m3 - We 
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t t 

(a) Total Cell Populations (b) Nutrient 

Figure 2: Results for the batch model, system ([1} with do = 0. We reculture 1% 
of each sub-population into new media at t = 70 hours. Functional forms are as 
in (|27| with parameters k = l/4hr, 7 = 4g CODB /m 3 , = 0.05, fx u — 0.005/hr, 
Y = 0.63g CODB /g CODS , C = gcoDB/5m 3 , k h = g COD B/(6hr m 3 ), and k g = 1/12 
hr. The dormancy domain is [0,1]. Time is measured in hours. 



require kh/C > l/24hr, so take kh — gcoDB/(6hr m 3 ). We set the dormancy 
domain to be so = and s* = 1, and assume 24 hr emergence at nutrient 
saturation so that k g = 1/12 hr. As we reculture or restore nutrients, we see 
a reawakening of the population and a growth spurt until the new nutrient 
is also depleted. Results are shown in Figure [2] for reculturing of 1% of each 
subpopulation onto new substrate at t = 70 hours. 

We conducted simulations with two species where 1% of each subpopula- 
tion was recultured into fresh media every 48, 72, and 168 hours. For shorter 
times of 4, 8, 16, and 24 (cases where 1% reculturing leads to extinction of 
all species), we used 100% reculturing. In all cases the fast-waker population 
(kg — 1/12 hr) outgrew the slow-waker population (k g = 1/24 hr). Although 
both active populations undergo oscillations, the fast-waker active population 
outgrows the slow-waker active population in each case. Moreover, the fast- 
waker total population dramatically outgrows the slow-waker total population 
in each case. (Results for 72-hour reculturing are shown in Figure |3]) This 
outcome is not surprising; conditions favor microbes that rapidly resuscitate. 

5.2 Chemostat Culture Dormancy As a Response to Nu- 
trient Depravation 

We use the functional forms and parameter values of Section [O] with do = k/2 
and Co(t) = 8 + 8cos(7rf/4). The active population, dormant population, and 
nutrient relax to a periodic oscillation. Results are shown in Figure 0J 



11 



Total Active Population 

Total Cell Population 

Total Active Population 

Total Cell Population 




B 5 
Z 

4 

3 
2 



..V 


\ 










































































































































I 




v 




I 




V 



150 200 



(a) Total Cell Populations 



(b) Nutrient 



Figure 3: Results for the batch model for two species with k g — 1/12 (blue) and 
kg = 1/24 (red). We reculture 1% of each subpopulation into new media every 
72 hours. All other parameters are as in Figure [2j 



For two competing species, the only difference being k g = 1/12 vs. k g — 
1/24, the results are shown in Figure [S] As in the batch case, the faster waker 
outcompetes the slower waker in the long run. This is true for a wide range of 
periods, C(f) = 8 + 8cos(27rtp) with p=0.5,l, 4, 12, 24, 48, 72 and 168 hours, 
verifying asymptotics predictions for short and long periods and extending to 
intermediate periods. Changes in period do not alter the numbers for a given 
subpopulation appreciably in magnitude, but rather change how they oscillate 
around some trajectory. The fast waker, as predicted, also outcompetes the 
slow waker in the case of a steady, rather than oscillating, nutrient source. 

Computations of competing species in a chemostat, where one species under- 
goes no dormancy (h = 0) and the other undergoes dormancy with parameter 
kg = 1/12, confirm the results of Section l3.2.1l (using C(t) = 8+8 cos(4-7rf) where 
t is measured in hours) and Section [3.2. 21 (using C(t) = 8 + 8cos(7rf/4) where t 
is measured in hours). A species without dormancy capability will outcompete 
an otherwise similar species which can go dormant, under both fast and slow 
oscillations in nutrient. 

5.3 Biofilm Dormancy as a Response to Nutrient Depra- 
vation 

In the biofilm model, (|27f| becomes 

f(c,u,d) = -7^) (" + / (e- s/k °+e- {1 - s)/k °)V ds^j . (28) 
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Figure 4: Results for the chemostat model, system ([T]) with do = k/2 and 
C(t) = 8 + 8cos(7rf/4) where t is measured in hours. Functional forms are as 
in (f27f with parameters k = l/4hr, 7 = 4g CODB /m 3 , = 0.05, (x u — 0.005/hr, 
Y = 0.63g CODB /g CODS , C = g C ODB/5m 3 , k h = g CO DB/(6hr m 3 ), and k g = 1/12 
hr. The dormancy domain is [0, 1]. Time is measured in hours. 



We use the functional forms and parameter values of Section 15.11 with the 
addition of p* = 10 4 g CODB /m 3 and D = 10- 4 m 2 /day [19]. 

Results shown in Figure |6] indicate that the slow waker has comparable 
total live biomass than the fast waker, and possibly more in lower regions, and 
palpably more dormant biomass, even though the faster wakers outnumber the 
slower wakers near the top of the biofilm. 



6 Computational Methodology 

As was done for senescence-structure in [SJ [T2] , we often handle a general physi- 
ologically structured system such as (jTa]) - tTTcl) more easily if it is transformed to 
an age-structured system, whether in the statement of the problem, or indirectly 
in the numerical method [5]. Let a > represent time a cell has spent dormant. 
We make a change of variables so that dormancy, s(a, t), is a separate function 
of age and time. We then get age-structured equations for v, 

d t v(a,t) + d a v(a,t) — —p, v (s,c)v(a,t) — Dv(a,t), 0<a<a*,t>0, (29a) 

v(0,t) = h(c)u(t), t>0, (29b) 

d t s(a,t) + d a s(a,t) = g(s(a,t),c(t)), 0<a<a*,t>a, (29c) 

s(0,t)=s o , (29d) 

where p, v (s,c) — p v {s,c) + d s g(s,c) and s(a*,t) = s* . The condition ( |27g[ ) 
means we need only consider t > a for the domain of s(a, t). 
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Figure 5: Results for the chemostat model for two species with k g = 1/12 (blue) 
and kg = 1/24 (red). All other parameters are as in Figure [4] In the long run, 
the slow waker is driven to extinction. 



For our choice, g(s, c) 



7g+C 



s = s 



for scalar c, we get 
Ig + c i T ) 



kgc{r) 



dr. 



(30) 



so that a* — > oo if c — > 0. Since functions with similar behavior to g are 
natural representations of the dormancy dynamics, we find that the original 
physiologically structured system is more tractable computationally than the 
equivalent age-structured system for most forms of g that interest us. 

To solve equations with more general physiological structure, we use an ex- 
tension of the natural-age-grid Galerkin methods developed for age- and space- 
structured systems in [21 [3] . These methods move the discretization nodes in 
age smoothly along characteristic lines. The solutions are approximated by 
piecewise polynomials, rather than moments as was done by de Roos [§]. Our 
extension of our methods to general physiological structure moves the discretiza- 
tion nodes in the physiological variable along characteristic curves, similar to a 
method of Sulsky |18] , but with the preservation of the property in our meth- 
ods that each time step need not result in a new discretization node in the 
physiological variable. This is essential when variation of spatial structure, or 
any other dynamics in the problem, occurs on a faster time scale than that of 
the physiological trait. Otherwise, the need to take lots of small time steps 
would induce many more physiological nodes than are necessary for accuracy, 
resulting in potentially great loss of efficiency from additional computation or 
interpolation onto a coarser grid. 

To motivate the integration in age and time, we ignore for the moment the 
discretization in space. We partition the domain [sq, s*] at each time by the set 
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Time=240 hours Time=240 hours 




cumulative volume fraction cumulative volume fraction 

(a) kg = 1/12 hr (b) k g = 1/24 hr 

Figure 6: Results for the biofilm model, system (|26|) . Functional forms are as 
in (p?T)) except that / given by ([28)1 . Parameters k — l/4hr, 7 = 4g CODB /m 3 , 
e h = 0.05, fj, u = 0.005/hr, Y = 0.63g CODB /g CODS , C = gcoDEs/Sm 3 , k h = 
g CODB /(6hr m 3 ). The two figures differ in the values of k g . The dormancy 
domain is [0, 1]. Time is measured in hours. The horizontal width of a color 
constitutes the volume fraction of cells of in the corresponding state. Green 
denotes active cells, blue denotes fully inert cells, and the "hot" black-red- 
yellow- white spectrum denotes dormancy from to 1. While fast reactivators 
are more prevalent near the top of the biofilm, slow reactivators have more 
biomass, particularly in dormant (hence resistant) form, in deeper biofilm layers. 



of nodes {si(t)}f =0 where So(t) = sn. If sjy(t) > s*, we simply ignore that node 
and the function value over it until needed. This is not an issue for our choices 
of g. We compute the solution at times tj and let Atj = tj+i — tj, s^j = Si(tj), 
and As,; j = s i+1 j — s,*j. For the last interval we use Asjv.j = max(s* — Sjv.jj 0)- 
Although we are not including space in this discussion for reasons of clarity, the 
presence of spatial structure in a problem will induce different time scales into 
a problem, making adaptivity and nonuniformity of time intervals an important 
property of any method used. 

For the computations in this paper, we use a piecewise constant approxi- 
mation space over the domain [so, s*]. Higher-order approximation spaces can 
be used, as was done in age in [2]. We define the projection into the space of 
piecewise constants over the partition of [sq, s*] by 

U(v(s,t 1 )) = \ As « Js *.i K ' 3> ' 1,3 - l+1 ' 3 ' (31 

v v ' 3 " \ 0, otherwise. v ; 

We make the approximation V^j IL(v(s,tj)) via the following algorithm. 
Let As max be the largest we want the first interval in s to be. For most time 
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steps we have Asqj < As max . In this case we set 

= s itj + Atjg(s it j,c(tj)), for i = 1,.. .,JV. (32a) 

We choose Atj such that Sij+i > s* for at most one i, so as to keep N fixed. 

Let the value Vij denote the density over [sij, Si+ij] for i = 1, . . . , N — 1. 
We use to denote the creation of newly dormant cells at sq. Then 

Vi,j+i = A Si,i V id , for t = 1, . . . , N - 1, (32b) 



As, 



3+1 



Fo J+ i = (AsojVoj + At 3 B(t 3 )), (32c) 

Because the applications in this paper provide for the first extension of the 
methods presented in [21 12], we have kept B{t) general in this part of the pre- 
sentation of the method. In our case we have B(tj) = g(so, c(t))v(s , t). Also, 
if g is independent of s, we have As.; j/Asj J+1 = 1 for i = 1, . . . , N — 1. 
If Asoj > As max , we introduce a new node and set 

Si+i,j+i = s itj + Atjg(sij,c(tj)), for i = 1,...,N- 1, (32e) 

Vi+ij+i = -^-V id , £ori = 0,...,JV-l, (32f) 

for the intermediate intervals, and set 

sij+i = At i ff(so,c(t i )), (32g) 



As 



3+1 



T/ A Sjv „ lj T/ jV _ lj + (Asjyj - At J - g (s*,c(t J -)))Viv, J - 

Kw.j+i = r , (32i) 

for the first and last intervals. 

The above calculations account for transport in the physiological variable, 
entry into dormancy, and exit from dormancy. Upwind differences approximate 
the advection terms in space. Center differences approximate the diffusion terms 
in space. Backward Euler formulae, embedded in step-doubling with local ex- 
trapolation, approximate the time derivatives. This creates a likely second-order 
correct time integration scheme [4]. 

7 Conclusions 

Modeling results suggest that spatial heterogeneity in biofilms can support a 
rich dormancy structure. For example, whereas dormant cells near the top of 
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a biofilm would need to be able to resuscitate quickly (small a*, large g) when 
environmental conditions improve in order to be competitive, dormant cells 
lower in the biofilm, where the slower waker has a defensive advantage over the 
fast waker due to a larger amount of dormant biomass without an appreciable 
difference in total live biomass, may be able to afford to be more cautious (large 
s* , small g). 

In contrast, dormancy-capable cells in well-mixed, planktonic systems (e.g. 
batch and chemostat cultures) appear to have less advantage over "regular" 
cells. In the absence of spatially structured populations, live biomass is max- 
imized by the fastest possible exit from dormancy. The lower limit of time 
to reawakening is governed by physiological, biochemical or other constraints 
within the cells, and hence dormancy mechanisms are constrained to easily re- 
versible mechanisms. As most lab populations are of the well-mixed batch or 
chemostat sort, and most natural populations are of the spatially-structured 
biofilm sort, this presents a possible drawback in use of typical laboratory sys- 
tems for characterization of natural ones. 

We remark that we have only considered here dormancy response in the 
context of resource depravation. Dormancy is also likely an effective defense 
strategy against antimicrobial agents - many antimicrobials are only effective 
against metabolically active targets. Thus the presence of antimicrobials rein- 
forces the utility of dormancy in biofilms and also may advantage dormancy- 
capable populations in well-mixed cultures. The nature of dormancy as defense 
could itself benefit from modeling studies. 

More generally, beyond dormancy specifically, recent studies suggest that 
phenotypic heterogeneity of many sorts is typical in spatially structured micro- 
bial populations such as biofilms [6]. Hence, methods of the sort presented here 
are likely to be useful and possibly necessary for modeling the function and 
ecology of spatially unmixed microbial populations of the sort that dominate 
the natural environment. 
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